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Abstract 

We use the method of solving the three-anyon problem developed in our earlier publica- 
tion to evaluate numerically the third virial coefficient of free anyons. In order to improve 
precision, we explicitly correct for truncation effects. The present calculation is about three 
orders of magnitude more precise than the previous Monte Carlo calculation and indicates 
the presence of a term asin^vrz/ with a very small coefficient a ~ —1.65 x 10^^. 
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1 Introduction 



The three-anyon problem, by now, has a rather long history (see and references therein). 
It is quite interesting as an example of a two-dimensional three-body problem which is 
apparently integrable, and besides that, it provides a means of determining the third virial 
coefficient of anyons. The virial expansion gives the pressure P in terms of the particle 
density n, 

Pp = n[l + A2i\^n)+A3{\\f + ...] (1) 

(/3 is the inverse temperature, A is the thermal wavelength) and shows the deviation of the 
equation of state from that of a classical ideal gas. Here A2, etc. are the (dimensionless) 
virial coefficients. For anyons, they depend on the statistics parameter u. 

The second virial coefficient of anyons, which can be found from the solution of the 
two-anyon problem, is known exactly 

Mi^) = 4 - ^^-^ , < z/ < 2 (2) 

{u = and 2 correspond to bosons, u = 1 to fermions). 

The problem of finding the third virial coefficient amounts, in one way or another, to 
solving the three-anyon problem and is therefore non-trivial. It is convenient to consider 
the difference 

AAsiu) = A,iu) - A^iO) , (3) 

where 

^3(0) = ^3(1) = ^ (4) 

is the third virial coefficient of bosons and fermions. On the basis of numerical calculations, 
it has been conjectured p, ^ that AA^lu) obeys a remarkably simple formula, 

AA^iu) = sin' Tiu . (5) 

It is consistent with the second-order perturbative result AA^lu) ^ z/^/12 P], with the 
exact supersymmetry property AA-^lu) = AA^ll — jy) ^ (which means that it is enough 
to consider u < 0.5), and with the numerical Monte Carlo results 0, ^ which are accurate 
to about 5% for u close to 0.5 . An independent check of (|^) by means of a numerical 
calculation of the energy levels has been carried out [0, |], confirming it for u < 0.25 with 
an accuracy of about 1% (the accuracy of the method being worse beyond this interval). 

This letter is devoted to a calculation of the third virial coefficient by finding the levels 
numerically using the method developed in our earlier paper [Q. With the inaccuracy in 
AA^iiy) estimated to be of the order of 0.005% for u close to 0.5, our results fit well to the 
modified formula ^ 

AAs{u) = — — - sin^ + a sin^ ttz/ (6) 

with a small coefficient a = — (1.652±0.012) x 10^^. This means that the task of calculating 
AA^lu) analytically appears more difficult than if (Rl) were exact. 
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2 The method 



Recall some basic equations concerning the third virial coefficient. It is most convenient 
to use the harmonic oscillator regularization, whence 



-21im 



where 



(7) 



(8) 



are the one- and two-particle partition functions, respectively, Z^ = Z^/Zi is the three- 
particle relative motion partition function, and x = h(3uo {uo is the harmonic frequency). 
The three-anyon spectrum consists of two types of states, interpolating linearly or nonlin- 
early between the bosonic and fermionic points, thus Z^ can be split into two parts. 
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^3 = ^3+ Z.^ 

where the contribution of the linear states is known exactly. 
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(9) 



(10) 



and the contribution of the nonlinear ones, Z^^, has to be evaluated. Since there is the 
so-called tower structure of the spectrum all the states coming in towers with the 

same angular momentum and with energy spacing 2hu}, it can be expressed as 

Z 



Z^' 
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where Z = Z{u,x) is the contribution of the bottom states (lowest states in the towers). 
Substituting (|])-(0) into (0) and then into (^ leads to the expression (cf. 0) 

Z(z/, x) — Z{0, x) 



-6 lim 



-2x 



We compute Z^u, x) according to its definition. 



Ziu, x) 



-xEsiiy) 



(12) 



(13) 



where s labels the nonlinear bottom states, and Es{i') are their dimensionless energies, 
obtained by setting huj = 1. These energies are found numerically making use of the 
method of the essence of which is as follows. The complex particle coordinates Zj = 



(X. 



+ iyj)/^/2, j = 1, 2, 3 are transformed into Z, u, v according to 





where rj = exp(27ri/3). The center-of-mass coordinate Z separates from the relative coor- 



dinates which are further rewritten as 



u 



rqe 



re 



i{</'-C/6) 



2(1 + g2 



2(1 + g2 



A common eigenstate of the relative Hamiltonian and angular momentum sat- 
isfying the cyclic interchange anyonic boundary condition, is ip{r,q,ip,^) = 
r^Pn{r'^) exp{—r'^/2)Q{q,(p,^), where P„ is a certain polynomial of degree n and the an- 
gular part is 

oo 

(](g,y.,0 = e'^^ E ^rngM^-'^"^^"'^'^^ (14) 

m=—oo 

with 

9M = + <ir"\F, (M±M _ I , - f ; 1 + bl ; -g^) , (15) 

L ( u'\ , L / u'\ 

?' = 3 mH , k = |-3mH 

2 V 2;' 2 ^ 2/ 

We have put uj and the particle mass equal to 1. The angular momentum is L = Lo + 3z/, 
where Lq must be an integer, v' = u + Lq mod 2, and the energy eigenvalue is 

E = ^i + 2n + 2. (16) 

The radial quantum number n is responsible for the tower structure, and the bottom states 
have n = 0. In the angular part of the wave function, the coefficients 7^ have to be chosen 
in such a way that it would satisfy the two-particle exchange anyonic boundary conditions, 
which have the form |^ 

00 00 

E 7m^7„^(l)e''"« = e'('^--'^'«) ^ 7,n^7^(l) e'^^ , (17) 

m=— 00 m=— 00 

00 00 

E 7^^7:n(l)e''"« = -e^('^--^'«) E l^Ci^)^-'"^^ (18) 



m=— 00 m,=— 00 



[ t/m(?) = dgm{q)/dq]. 

To use these boundary conditions, introduce a linear operator A such that 

[^0](e)=e'(— '«)</,(27r-0. (19) 

We take the function (f){^) to be defined for ^ G [0, 27r] and we represent it by the vector of 
its Fourier components {0m}, so that 0(^) = J2m=-oo 4>mG^"^^- The matrix elements of A in 
this representation are 

sinvri/ 

= ; ; jr ■ (20) 

Next introduce two diagonal matrices 

Gmn 9m(X) ^mn j ^mn 9nSX) ^mn ■ (21) 

The two boundary conditions ([T7|)-(p!8|) then take the form 

(/ - A)G^ = , (/ + A)G'^ = , (22) 

where 7 is the vector {7™,} of the unknown coefficients in the wave function ([T^). Now, 
since A is a real symmetric matrix and = J, the vectors on the left-hand sides of the 



two equations in (p2D are orthogonal. This makes it possible to replace the two equations 
by one, 

[G + G' + A{G' - G)]7 = . (23) 

A nontrivial solution for 7 exists if and only if the determinant, which depends on /i through 
(pip and (|T3p, vanishes. 

Thus, in order to find the energy eigenvalues, one fixes chooses a sector of bosonic 
angular momentum Lq and scans an interval of /i in order to find the zeros /is(z^) of the 
determinant. After repeating the procedure for different Lq and dropping the known linear 
states, one gets the eigenvalues Egiy) = Hsiy) + 2, to be substituted into (plsl). A good 
consistency check of our numerical results is that the number of states found agrees with 
the exactly known multiplicities [llO . 



3 Truncation effects 

There are two sources of error in our calculation even at infinite numerical precision. The 
first one is only a finite number of Fourier coefficients being considered, the second one is 
only states with energies less than a certain finite £ being included ("energy cutoff"). The 
two effects are independent, and in order to improve precision we will explicitly take them 
into account. 

We truncate the Fourier series at \m\ = M, i.e., consider {2Af + 1) Fourier compo- 
nents {7__yv', 7-A^+i, • • • ,7Ar}. To maintain consistency, one has then to restrict oneself to 
a discrete set {^__^/-, • • • ,^m} instead of a continuous variable ^, so that any set of 

function values 7(^n) can be represented by the Fourier components exactly. Choosing 



= + 27rr7,/(2A/' + 1) and demanding (|T9[) to hold at the points C,n, one gets for the 
matrix elements of A 

A = ^^I^ (24) 

™" (2A/'+l)sin[7r(m + n + i/')/(2Ar+l)] ' ^ ^ 

with the correct Af ^ 00 limit . 

By analyzing the data Es{i',M) obtained at different A/", we find that they are rather 
well described by an empirical formula 

E.(.,Ar) = E.(.) + -A + _lL_ + ... , (25) 

which allows to extrapolate to M = 00. 

Now assume that the levels are found exactly and consider the effect of the energy 
cutoff. In practice one can only find the quantity 

Z{v,x,£)= Y: e--""^^^^ (26) 

Es{u)<£ 

for a finite £, so that ([I^ ) should be rewritten as 

AA^iu) = lim lim AAs{iy,x,S) , (27) 

x~*0 £—>€<> 



where 

Z{u, X, 8) — 2{0, X, 8) 
1 -e 



AA3(., X, 8) = -6 ^^'^'^'r^:^r^'^^ . (2^ 



A 



Note that the order of hmits cannot be interchanged. Of course, Z{0, x) is known exactly, 
but truncating it at the same £ removes a spurious divergence. 
We now argue that 

AAs{u, X, £) = AA3(z/, x) + e-^V(i^, S) , (29) 

where /(z/, x, grows (or diminishes) with x not faster than e^^ with \k\ ~ 1 and has a 
weak (polynomial type) dependence on £. Indeed, the number of nonlinear bottom states 
grows with energy E like i?^/9 [0], and, respectively, 



9(1 -e^')(l -e^"^)^ 

(up to terms proportional to £ mod 3). Since when v changes, energies change by values 
of the order of unity but the number of states is conserved, Zi^u, x, £) should still have the 
form Zi^u, x) + e~^'^JF(z/, x, £) with JF changing not faster than £^ and not faster than . 
This leads to the statement expressed by eq. (|29|). In fact, will be singular at x — for 
a finite £, but all the singularities cancel out when doing the subtraction in (p8|). 
On the other hand, at small x there has to be an expansion 

Ay43(i/, x) = lim AA^i^u, x, £) = ao + a2X^ + a^x^ + . . . . (30) 

For a discussion of why we include only even powers, see the next section. Substituting 
(pop into (p9D and making a power series expansion for f{i',x,£), one has that 

AAs{i', X, £) = ao + a2X^ + a^x'^ + . . . + e"^^(6o + hx + h2X^ + ...). (31) 

Remember that for any finite £, AA^{v,x^£) is finite at x = 0; this is in contrast 
to ref. where it was treated as being infinite. In fact, expanding the numerator and 
denominator of (|28|) in powers of x, one gets 

AA,[v,x,£)- 6^c»^^^_-^)„+i^2x)Vn! 

= 35i + (^35i - ^^2) x+[s,- ^2 + iSs) x^ + . . . , (32) 

where 

Sn^s^{u,£)= y: [E:iu) - e:{q)] (33) 

Esiu)<£ 

is finite for any n. The ansatz (|3TD must be made consistent with (0). Matching coeffi- 
cients at equal powers of x after expanding the exponential in (^) leads to the following 
expressions for the coefficients bn- 

bo = 35*1 — , 

3 

6i = (3 + 3£)Si - -S2 - £ao , 

1 + 36 + -— b^l - 77 + P2 + 7:03 - «2 , 
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This leaves only oq, 02, • • • as unknown coefficients. 

We now retain only a finite number of coefficients oq, 02, • • • , o,2n and bQ,bi, . . . , in 
(pT|). The truncated equation can be rewritten as 



N 



n=0 



(34) 



with all the unknown coefficients on the right-hand side and all the terms depending on 
the energies Es{v) on the left-hand side. The quantity on the left-hand side is 



2N 



n=0 



with 



Co = "iS 



1 ! 



ci = (3 + 3^)5i--^2, 



(35) 



and the functions on the right-hand side depend only on x and S, 

2N-2n 
m=0 



f2nix,S) = X 



2n 



(36) 



From eq. (Q) with AA'^{iy, x, £) calculated numerically, one obtains the coefficients a2n 
through a least square fit over a certain x interval, and by virtue of (|27| ) and (^) one has 
AA^lu) = aQ. Note that at a; = 0, each of the functions f2n{x,£) vanishes together with 
its ffist 2N — 1 derivatives with respect to x. Indeed, by the subtraction in ( pSj ) one has 
already taken into account the small x behavior of AA^lujXjS), so that AA^{i>,x,S) has 
to vanish at x = up to order 2N — 1. But fixing a positive x and taking £^ — 00, as in 



results in f2n{x, S) 



X 



2n 



meaning that AA3 (z/, x, 8) tends to the 8 = 00 form 



4 Results and discussion 

We take 6 different values Af = 10, 20, 40, 80, 160, 320 and the energy cutoff up to £^ = 35. 
The Af convergence is better for u closer to 1 — which is natural, because then the wave 
functions f2(l,0, ^) are less singular at ^ = 0,27r — and therefore we do all the calculations 
on the (semion-fermion) interval u G [0.5,1]; that is, the values listed for u are actually 
obtained at 1 — z/. Taking six different values of Af allows to include five correction terms 
in dH). 

The procedure of extracting AA^lu) is illustrated in Fig. 1, for the case of semions, 
u = 0.5. (In all the figures, the data for AA3 are multiplied by 10^, for convenience.) 
The X expansion was truncated at 2N = 8. As £^ is increased, AA^{i',x,8) approaches 
the 8 = 00 curve at any finite x, but always vanishes at x = 0. We emphasize that, 
unlike in ref. |]^, we do not need an x cutoff on the left; the fit is performed on an interval 
X G [0,Xmax], and the points at small x, where the initial data have nothing in common 
with the 8 = 00 curve, are fully accounted for. 



We take different £ and Xmax and observe how the extracted ^A^^v) depends on those. 
Ideally, there should be no dependence at all, thus the dependence left over provides an 
estimate of the error of the calculation. Figure 2 shows this dependence for semions. For 
a^max too small, the results of the fitting become senseless, since one tries to extract too 
much information from too narrow a range of data, while for Xmax too large, the effect of 
truncation of the x expansion becomes significant. It is seen, however, that there exists 
an interval of Xmax where /S.A'i{v) is practically constant. Shown in the figure is an "error 
bar" that we take for the result. Increasing or decreasing 2N leads, in accordance with the 
aforesaid, to the "good" interval of Xmax moving to the right or to the left, respectively, 
but the result does not change significantly for 2N changing down to 4 or up to 10. On 
the other hand, the £ and a;max dependence becomes considerably stronger when the data 
for M = 160 and/or 320 are not used. Thus, the main source of error appears to be in the 
inaccuracy of individual levels. 

We want to emphasize that our error analysis by varying the parameters Af, £ and 
Xmax is at the same time a justification of the methods we use for obtaining the limits as 
A/" — > oo, £^ — i> oo and x — 0. We argue that the result for AA^i^u) is correct to the degree 
to which it is independent of £, x^aax, and the maximal value of A/". As for the limit x — 
in eq. (P7|), if eq. (^TJ) did not represent the true x dependence, we would expect the fitted 
values of the coefficients a2„ in eq. (|3^) to depend on the interval < a; < Xmax where the 
fit is done. Since there exists a finite range of Xmax where they depend very little on Xmax, 
we believe that these are close to the true values. 

The omission of odd powers of x in eq. (^) (like in ref. |^) is only empirically justified, 
apart from a general argument |Tl| to the effect that the linear term must be absent, as this 
expansion results from an h expansion of the partition functions (we thank Avinash Khare 
for pointing this out). As a second check, we tried to include the terms x^,x^, .... This 
changed the values of AA^^u) by amounts typically about twice our quoted errors, while 
the coefficients of these terms were one to two orders smaller than those of the adjacent 
even order terms, and were minimal for values of x^ax where AA^lu) depends most weakly 
on Xmax- This is a clear indication that these terms are really absent, so that the fit where 
they are omitted should be more precise. 

Table 1 and Fig. 3 contains the main result of this paper. Shown there is the difference 
AA^i^u) — (1/1277^) sin^ ttu together with the best fit curve of the form a sin^ vrz/. The result 
is 

a = -(1.652 ±0.012) X 10"^ = . (37) 

^ ' (621±5)7r4 ^ ^ 

The value of per degree of freedom is 1.0. If we exclude the three points closest to u = 0, 
it drops to 0.45, which we regard as a strong indication that the errors are overestimated 
for u close to 0.5. This is in contrast to ref. where the results for u close to 0.5 were 
the least accurate. We have included the value z/ = (3 — v^)/2 ~ 0.38, the golden ratio, 
as an example of an irrational number which is badly approximated by rational numbers. 
This point falls on the same smooth curve as all the points at rational values of u. 





Our result 


Eq. (§ 


0.05 


2.063(2) 


2.066 


0.10 


8.057(3) 


8.063 


0.15 


17.391(3) 


17.403 


0.20 


29.149(2) 


29.171 


0.25 


42.174(2) 


42.217 


0.30 


55.192(2) 


55.263 


0.35 


66.927(3) 


67.032 


0.38 


73.222(2) 


73.347 


0.40 


76.237(2) 


76.372 


0.45 


82.212(2) 


82.368 


0.50 


84.270(2) 


84.434 



Table 1. The values of WAA^{v). 



Shown in Fig. 4 for comparison are the Monte Carlo data obtained by the method of 
ref. from 330 million paths at a; = 0.25 and from 131 million paths at a; = 0.35 0. 
By construction, those data have a real part (drawn by solid line) and an imaginary part 
(dashed line). The imaginary part should vanish, hence its deviation from zero gives an 
estimate of the statistical error in the real part. The leading x dependent term in the latter 
is expected to be proportional to . Since 0.35^/0.25^ ~ 2, we see from the plots that 
the result of an extrapolation to x = would not be significantly different from zero. One 
sees that our present calculation is about three orders of magnitude more accurate than 
the Monte Carlo calculation. The sin^ tiv term that we find could not have been seen in 
that calculation, neither in the one of ref. [Q, because of insufficient precision. 

As originally noted in ref. 0, periodicity, supersymmetry and analyticity imply that 
the third virial coefficient can be represented as a rapidly converging Fourier series in only 
even powers of sinvrz/. The coefficient at the fourth power calculated here is about 500 
times smaller than the one at the second power. A way to find it exactly would be to do 
fourth-order perturbation theory in the manner of ref. 0. Our calculation is not precise 
enough to show whether there are higher powers of sinvrz/, in fact the Fourier series might 
well be infinite. 

The consistency of our present results with all previously known results is another good 
check that they are reliable. In particular, that we get a sin^ vrz/ term with the coefficient 
known from perturbation theory is a nontrivial result. That the deviation from this sin^ vrz/ 
term has the form sin'^ vrz/ is equally nontrivial. If we try to fit it as sin^ vrz/, for example, 
we get xVDOF = 80. 

It is a remarkable result that the conjecture (|^) is valid within about 0.2%, even though 
it is not exact. This is likely to mean that there exist some, as yet undiscovered, approxi- 
mate analytic formulas for the nonlinear energy levels, that yield formula (|^) and are valid 
within a fraction of a per cent. It would be rather interesting to understand what those 
formulas are. 

We thank the Centre for Advanced Study for kind hospitality and financial support. 
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Figure captions 



Fig. 1. Extracting AA-^lu) from the data for S = 35, at u = 0.5 (semions). 
Triangles: 10'' AA3(z/, x, eq. (p^), versus x. 
Diamonds: 10^ AA^ {u, x, S), eq. ®. 
Solid line: the fit by eq. (0). 

Dashed line: the S = oo curve, eq. (|30D . The a; = point of the latter is AA^lu). 

Fig. 2. The Xmax dependence of the extracted 10'^ AA^i^u) aX v = 0.5. Curves for £ = 25,30 
and 35 are plotted. Also shown is our final result with an error bar. Note that 
10Vl27r2 = 84.43 is well outside the plot. 

Fig. 3. The difference 10'^ [AA3(z/) — (l/127r^) sin^ vri/] as a function of u. The points with 
error bars are our present results. The curve is a least square fit. 

Fig. 4. The same quantity as in Fig. 3, from the Monte Carlo simulation. 
Solid lines: real part. 

Dashed lines: imaginary part (illustrates the size of statistical errors). 

The upper part is for x = 0.25 (330 million paths). The lower part is for x = 0.35 

(131 million paths). 
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